#include <mln/all.hh>
#include <list>
#include <utility>
#include <cmath>

typedef enum
  {
    DIR_RIGHT = 0,
    DIR_RIGHT_UP = 1,
    DIR_UP = 2,
    DIR_LEFT_UP = 3,
    DIR_LEFT = 4,
    DIR_LEFT_DOWN = 5,
    DIR_DOWN = 6,
    DIR_RIGHT_DOWN = 7,
    DIR_NONE = 8
  } e_dir;

template <typename T>
static inline
T*
next_neigbh(int ncols,
	    T* p,
	    e_dir *dir)
{
  switch (*dir)
    {
    case DIR_RIGHT:
      return p + 1;
    case DIR_RIGHT_UP:
      return p + 1 - ncols;
    case DIR_UP:
      return p - ncols;
    case DIR_LEFT_UP:
      return p - 1 - ncols;
    case DIR_LEFT:
      return p - 1;
    case DIR_LEFT_DOWN:
      return p - 1 + ncols;
    case DIR_DOWN:
      return p + ncols;
    case DIR_RIGHT_DOWN:
      return p + 1 + ncols;
    default:
      return p;
    }
}

static inline
void
dec_dir(e_dir* d)
{
  *d = (e_dir) (((unsigned)*d + 7) & 7);
}

static inline
void
inc_dir(e_dir* d)
{
  *d = (e_dir) (((unsigned)*d + 9) & 7);
}

static inline
void
inverse_dir(e_dir* d)
{
  *d = (e_dir) (((unsigned)*d + 4) & 7);
}

template <typename T>
static inline
e_dir
find_first_dir(unsigned ncols,
	       T* p)
{
  e_dir dir = DIR_NONE;

  T* n;
  do
    {
      dec_dir(&dir);
      n = next_neigbh(ncols, p, &dir);
    }
  while (*n != mln::literal::zero);

  while (*n == mln::literal::zero)
    {
      inc_dir(&dir);
      n = next_neigbh(ncols, p, &dir);
    }

  return (dir);
}


template <typename T>
static inline
T*
contour_next(unsigned ncols,
	     T* origin,
	     e_dir* dir)
{
  T* neigbh;

  do
    {
      dec_dir(dir);
      neigbh = next_neigbh(ncols, origin, dir);
    }
  while (*neigbh == mln::literal::zero);

  inverse_dir(dir);

  return neigbh;
}

template <typename T>
static inline
unsigned
hook_to_index(const T* hook,
	      const T* first)
{
  return ((unsigned) (hook - first)) / sizeof (T);
}

template <typename I>
static
std::list<mln::point2d>*
code_base(const I& ima,
	  typename I::value* origin)
{
  typename I::value* current = origin;
  typename I::value* old = origin;
  typename I::value* node = origin;
  const typename I::value* first = ima.buffer();
  unsigned ncols = 2 * ima.border() + ima.ncols();
  e_dir dir_origin = find_first_dir(ncols, origin);
  e_dir dir = dir_origin;
  e_dir old_dir = dir;

  std::list<mln::point2d>* l = new std::list<mln::point2d>();

  do
    {
      if (dir != old_dir)
	{
	  l->push_back(ima.point_at_index(hook_to_index(node, first)));
	  node = old;
	}

      old = current;
      old_dir = dir;
      current = contour_next(ncols, current, &dir);
    }
  while (current != origin);

  l->push_front(ima.point_at_index(hook_to_index(node, first)));

  return l;
}

template <typename I>
std::list<mln::point2d>*
code_base_initialize(I& ima)
{
  unsigned i = 0;
  typename I::value* start = ima.buffer();
  unsigned ncols = ima.ncols() + 2 * ima.border();

  start += ima.border() * (ncols + 1);

  i = 0;

  start += ima.ncols() / 2;

  while ((i < ima.nrows()) && (*start == mln::literal::zero))
    {
      start += ncols;
      ++i;
    }

  return code_base(ima, start);
}

template <typename I>
static
void
draw_polygon(I& ima,
	     std::list<mln::point2d>& l)
{
  if (l.size() == 0)
    return;

  mln::point2d prec = l.back();

  for (std::list<mln::point2d>::const_iterator i = l.begin();
       i != l.end(); ++i)
    {
      mln::draw::line(ima,
		      *i,
		      prec,
		      mln::literal::red);
      prec = *i;
    }
}

int main (int argc, char** argv)
{
  mln::image2d<bool> ima;

  if (argc != 2)
    return 1;

  mln::io::pbm::load (ima, argv[1]);

  std::list<mln::point2d>* contour = NULL;

  mln::util::timer t;

  t.start();


  contour = code_base_initialize(ima);
  t.stop();

  std::cout << t << std::endl;

  mln::value::rgb8 v;
  mln::image2d<mln::value::rgb8> res = mln::data::convert(v, ima);

  if (contour)
    {
      draw_polygon(res, *contour);  
      delete contour;
    }

  mln::io::ppm::save(res, std::string(argv[1]) + ".out");

  return 0;
}
